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1.  INTRODUCTION 


Though  assumed  to  be  time  homogenous  in  most  studies  of  queueing 
models,  arrival  processes  are  seldom  such  when  real  systems  are  exam- 
ined. In  addition  to  the  difficulties  involved  in  the  mathematical 
analysis  of  queueing  systems  having  non -homogenous  arrival  processes, 
it  appears  that  variations  in  these  processes  are  generally  slow,  rela- 
tively to  the  queueing  process  itself,  thus  time  intervals,  in  which 
arrival  rate  may  be  considered  constant,  are  long  enough  to  allow  steady 
state  analysis.  However,  often  this  is  not  the  case,  and  non-homogeneity 
of  the  arrival  process  must  be  assumed.  Queueing  models  with  non-homo- 
genous  input  have  attracted  the  attention  of  queueing  theorists  starting 
in  the  mid-fifties  with  the  works  of  Clarke  [2]  and  Luchak  [8].  Since 
then  there  is  a thin  flow  of  studies  of  such  models,  in  the  evergrowing 
stream  of  literature  of  queueing  theory  (see  references  [1]  to  [10]). 

Moore  [9]  divides  the  approaches  to  "non-stationary"  queueing  sys- 
tems into  two  primary  categories:  those  considering  the  parameters  of 

the  arrival  (and  service)  process  (processes)  to  be  continuous  functions 
of  time,  and  those  assuming  parameters  to  be  step  functions.  Yet,  there 
are  still  other  possibilities.  Many  systems  are  characterized  by  the 
random  nature  of  the  variations  of  the  basic  processes.  Hasofer  [3] 
describes  a situation  where  arrivals  come  from  a large  number  of  N of 
"Independent  sources",  each  having  a Poisson  output  with  parameters  \, 
hence  arrivals  are  of  Poisson  type  of  parameter  N\.  If  the  number  of 
sources  is  assumed  to  be  a function  of  time  N(t),  the  arrival  process 
will  be  a non-homogenous  Poisson  process  with  time  dependent  parameter 


\N(t).  Notice  that  in  this  presentation,  N(t)  is  assumed  to  be  deter- 
ministic. 

In  many  real  life  systems,  however,  the  variations  in  the  number 
of  sources  are  due  to  stochastic  phenomenae,  thus  modeling  of  such 
systems  calls  for  the  replacement  of  the  (deterministic)  function  N(t) 
by  a stochastic  process  ^(t),  t > 0].  Obvious  candidates  are  those 
describing  the  behavior  of  populations  (Poisson  process,  birth  and  death 
process,  etc.).  Restricting  ourselves  to  Markov  processes,  we  represent 
the  arrival  rates  by  a Markov  jump  process  fA(t),  t > O},  on  the  param- 
eters states  set  with  transitions  ^2  input 

parameter  process  is  discussed  in  Section  3* 

It  should  be  noted  that  step  functions  are  often  used  to  approximate 
continuous  functions,  with  proper  selection  of  jump  points  and  values. 
Similarly,  certain  (deterministic)  arrival  rate  functions  may  be  approxi- 
mated by  stochastic  processes  of  the  type  considered  in  the  present 
study,  with  proper  selection  of  the  parameters  state  set  and  transition 
rates. 

We  consider  here  Markovian  queueing  models.  For  reasons  of  simpli- 
city we  concentrate  on  the  M/M/1  model.  The  results  may  be  however, 
easily  extended  to  more  general  classes  of  Markovian  models,  defining 
birth  and  death  processes  with  stochastically  dependent  birth  intensities 
and  state  dependent  death  intensities.  For  example,  an  M/M/r  model 
will  be  specified  by  the  arrival  rate  process  [A(t),  t > O)  and  service 
rate  function  |i(x)  = min{nx/r,|a],  when  x is  the  current  queue  size. 

Queueing  systems  with  non-homogenous  input,  may  never  attain  sta- 
tionarity,  thus  steady  state  solutions  cannot  be  used.  On  the  other  hand, 
time  dependent  solutions  may  not  yield  measures  of  performance  which  are 
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essential  if  solutions  are  to  be  used  for  decision  making.  In  such 
cases,  either  expected  accumulated  costs  (for  finite  horizon  only), 
or  expected  accumulated  discounted  costs  should  be  considered.  Accu- 
mulated costs,  associated  with  decision  processes  will  be  discussed  in 
a following  communication.  In  the  present  paper  we  consider  discounted 
costs. 

Some  of  the  results  are  based  on  previous  results  by  Yadin  and 
Zacks  [11]  and  [12],  and  others,  the  summary  of  which  is  given  in  Sec- 
tion 2. 

The  main  results  are  presented  In  Section  k. 

2.  PRELIMINARIES 

A wide  class  of  Markovian  queueing  models  possess  the  following 
properties : 

(a)  input  and  output  processes  are  both  Markovian; 

(b)  actual  service  rate  (output  intensity,  death  intensity) 
is  uniquely  determined  by  the  current  queue  size.  It 
cannot  exceed  a given  bound  called  the  service  capacity. 

Considering  discounting,  a model  is  specified  by: 

(a)  a set  of  parameters  j\,n,p],  where  K,  \ > 0 is  the  arrival 

rate;  )i,  > 0,  the  service  capacity;  and  p,  p > 0 the 

discounting  rate. 

(b)  an  output  function  ii(x),  0 < n(x)  < n,  representing  output 
Intensities  as  a function  of  the  queue  size  x,  x = 0,1,...  . 

For  example,  the  output  function  of  the  M/M/r  model  is 
(2.1)  n(x)  = min{^lx/r  , n]  . 

The  original  set  of  parameters  can  be  replaced  by  where 
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(2.2)  § = + P ” + p)^  ” 4\n  ) ]f 

ijt  = + H + p+  ((>^  + H + p)^  - 4\|i)^^^]  . 

It  Is  well  known  that 

(2.3)  H = \ I ii 

p = \(l-§)  (i|r-l)  . 

Consider  a queueing  process  fX(t),  t > O},  X(t)  designating  the 
queue  size  at  epoch  t,  with  a transition  function 

(2.4)  P = y|X(0)  = x}  , x,y  = 0,1,...,  , 

depending  upon  {X,n,p},  (or  alteratively  {X,§,'lf})  and  n(x).  We 
assume  that  there  exists  a number  Xq  such  that  |^(x)  = n for  all 

X > Xq. 

For  any  P,  3 > 0,  let  0^  denote  an  independent,  exponentially 
distributed  random  variable  with  mean  1/p. 

Let 

(2.5)  m(x)  = f S m(y)  p_(t)  pe“^*=  dt 

0 

= E[m(X(e^)|X(0)  = x]  , 

be  defined  for  any  function  m(x),  x = 0,1,...  such  that  the  conver- 
gence condition 

(2.6)  11m  r*  ®(x)  = 0 

x-»0 

is  satisfied.  m(x)  defined  in  (2.5)  Is  known  as  the  3-potential  of 
the  charge  3o(x) . 
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It  is  easily  seen  that  is  a linear  operator,  satisfying 


n n 

(2.7)  S a m (x)]  = a + S a m (x)  , 

" i=l  ^ ^ i=l  ^ ’• 

for  any  set  of  coefficients  (aQ,...,a^]  and  charges  fPnij^(x), . . . ,3m^(x)'|, 

O 

such  that  Umj(x),  j = l,...,n,  exist. 


(2.8) 


G^(z,p)  = j 


FI  < ^ , 


be  the  Laplace  transform  of  the  generating  function  of  the  transition 
distribution  queueing  process  fX(t),  t > Ol,  evalu- 
ated at  p.  It  can  be  shown  that  for  | z|  < (^  be  a function  of  the 

discounting  rate  p)  the  convergences  conditions  are  satisfied.  The 
transform  G^(z,p),  which  is  given  for  a variety  of  models,  may  conve- 
niently yield  U^m(x)  for  a wide  class  of  charges,  in  the  same  fashion 
as  generating  functions  generate  probabilities  and  moments.  Indeed,  if 
(x)^  = x(x-l) • • • (x-r+1)  then 


(2,9) 


(2.10) 


-^(x)  z*  G (z,p)  , I z;|  < if  , 

P dz 


“ ^ G (z,p)|  ^ 

p r,x  r x'  I z=0 


where  6 =l,r=y:6  =0  otherwise.  In  particular, 

r,y  ' r,y  ^ ’ 


(2.11) 


= :4  V"'p>iz=i  • 


Clearly,  computations  based  on  the  last  results  are  practical  for 
small  Integers  r only.  Combining  these  results  with  (2.7)  we  obtain 
an  efficient  tool  for  the  evaluation  of  a m(x)  for  various  charges. 
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Further  analysis  requires  the  consideration  of  specific  models. 
For  reasons  of  simplicity  we  shall  concentrate  on  the  M/M/1  model. 
Extensions  to  other  models  are  immediate. 

It  is  well  known  that  for  the  M/M/1  model 


(2.12) 
particularly 


1 /I  ^■'■1  /I  N + 1 

uP  i = G (z,p)  = 

P ""  X(l-5)(z-5)(i(;-z) 


z|  < ♦ . 


(2.13) 


p 


g-g)  X +1  X 


It  follows  from  (2,7)>  (2.11)  and  (2.12)  that 


(2.14)  uP[s  + qx]  = - 

P 

as  shown  by  Yadin  and  Zacks  [11]. 

We  are  considering  in  the  present  study  the  conditional  expecta- 
tions of  accumulated  discounted  costs.  Let  M(x)  denote  such  expecta- 
tion evaluated  at  t,  t > 0>  given  X(t)  = x.  Because  of  the  Markovian 
nature  of  the  model  this  expectation  is  independent  of  t.  Assuming 
that  the  costs  are  generated  by  a total  charge  of  pm(x)  per  time  unit, 
we  have 


(2.15) 


M(x)  = uP  m(x)  . 


The  charge  pm(x),  consists  of  holding  costs  of  x customers  and  other 
possible  costs.  The  costs  function  m(x)  represents  the  accumulated 
discounted  costs  associated  with  the  residence  of  x customers  for  an 
Infinitely  long  period. 

Obviously,  we  shall  limit  ourselves  to  costs  functions  such  that 
M(x)  is  finite  for  any  finite  x.  Furthermore,  we  shall  assume  that 
m(x)  is  such  that  the  convergence  condition  (2.6)  is  satisfied.  Linear 
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1 X 

functions,  polynomials  and  exponential  functions  — z such  that 

P 

[z]  < (|f,  are  examples  of  such  m(x)  functions. 

The  process  {X(t),  t > 0}  is  a birth  and  death  process  with  birth 
intensity  \ and  death  intensities  n(x).  Transition  probabilities  are 
+ n(x)]  ^ and  p(x)[\  + n(x)]  respectively.  Let  0 denote  the 

random  time  of  the  first  event  (either  birth  or  death).  0 has  an 
exponential  distribution  with  parameter  \ + n(x).  It  follows  that 


(2.16) 


M(x)  = E j J pm(x)  e ^ dt  + 


+ exp(-p0) 


|x.+ia(x) 


M(x+1)  + 


^ M(x-U 


Integration  of  (2,l6)  yields  the  following  functional  equation  for  M(x) : 


(2.17) 


M(x)  = 


\ M(x+1)  + u(x)  M(x-l)  + 
\ + n(x)  + p 


The  unique  solution  of  the  functional  equation  (2.17)  which  satis- 
fies the  convergence  solution  is  the  one  presented  at  (2.15).  Suppose 
there  exist  two  distinct  solutions,  M(x)  = U^m(x)  and  M(x)  + d(x). 
Clearly,  the  difference  d(x)  satisfies  the  equation. 


(2.18) 


d(x)  = 


LX+j.)  + u(x;  d(x-J 
\ + n(x)  + p 


, for  all  X. 


Hence,  for  all  x > x^  (for  which  u(x)  = |j) 


(2.19) 


d(x)  = a I + b i|f 


In  order  to  satisfy  (2.6)  we  must  set  b = 0. 

x-x 

Substitution  for  M*(x)  and  M*(x)  + a | into  (2.17)  (for 

X > 0)  yields 


(2.20) 


a[\5(l-|)  - n(l-5)  + p]  = ap(l-P)  = 0 . 


Since  p > 0 and  ? < 1 we  must  have  a = 0.  Thus  d(x)  = 0 for  all 


X > Xq  and  therefore  for  all  x. 

M(x)  can  be  approximated  by  a sequence 


such  that 


(2.21) 

It  can  be  shown  that 
of  (2.17),  for  any  M^^^(x) 
not  violated. 


\M^*^'^^(x+l)  + u(x)  M^^~^^(x-1)  + pm(x)  . . 

\ + H(x)  + p ' - 

(M  (x) } uniformly  converges  to  the  solution 
such  that  the  convergence  conditions  are 


3.  THE  MAIN  RESULTS 

Let  aQ,aj^,...  and  be  two  sequences  of  non-negative 

parameteis,  and  let 


(3.1) 


■^0  = ° 

" Vi  ’^k-i  ’ ^ 


where  Tq,T^,...  are  mutually  independent,  exponentially  distributed 
random  variables  such  that  E[T^]  = ctj  , j = 0,1,...  , 

Define  a Markov  process  (A(f)  > f ^ 0],  such  that 


(3.2)  A(f)  - ^ ^ ^k+1  ’ k — 0,1,2,..,  . 

By  setting  Q = 0 for  some  n we  mean  that  A(f)  = ^ < t < «, 

n n n 

in  other  words,  the  last  change  in  A(t)  Is  at 

Consider  a non-homogcnous  Poisson  process  TA(t)  , t > 0)  with 
intensity  A(t)  at  t,  t > 0.  The  objective  of  the  present  work  is 
the  evaluation  of  the  expectations  of  the  accumulated  discounted  costs 
of  queueing  systems  of  the  class  presented  in  Section  2,  with  the  non- 
homogenous  arrival  process  [A(t)  , t > 0"}.  Note  that  the  arrival  process 
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is  completely  defined  by  a = faj^}  and  \ = Hence,  the  model  is 

specified  by  the  set  of  parameters;  by  the  service  function 

p(x)  and  by  the  costs  function  m(x). 

Denote  by  conditional  expectation  of  the  accumulated 

discounted  costs  evaluated  at  some  epoch  t,  for  infinite  period,  given 
Tj^  < t < and  X(t)  = X. 


(3.3) 


00  00 

Mj^(x)  = J S p m(y)  P [X(t)  = ylx(t) 
, y=0 

. e-P^^-'^^dT  = 


= E[m(X(t  +0p))lX(t)  = X,  Tj^  < t < . 

Arguments,  similar  to  those  leading  to  (2.16)  and  (2.17)  yield 


(3.^) 


0 


-pt 


Mj^(x)  = E < J p m(x)  e dt  + 


+ e-P® 


+ M,(x) 


\(x+l) 


a. 


^ Xj^  + + u(x)  \^+  a^  + ^i(x) 


where  9 has  an  exponential  distribution  with  parameter  ^ l^(x). 

This  expression  is  reduced  then  to 


(3.5) 

Following  the 
solution  of  (3.5) 


\ ^(^■•■1)  l^(^)  Mj^(x-l)  + M^^j^(x)  + pm(x) 

“ Xj^  + (i(x)  + + p 

results  of  Section  2,  we  observe  that  the  unique 
satisfying  the  convergence  conditions,  is 


+P 

Mj^(x)  = U f 


(3.6) 


where 


(3.7) 


mj^(x)  = 


ak  + P 


Given  the  costs  function  m(x)  and  some  n > 0,  ^(^) 

may  be  recursively  determined  for  all  0 < k < n,  either  numerically 
(applying  approximating  sequences  generated  by  equations  similar  to 
(2.21)  and  based  on  (3-5))  or  analytically. 

Clearly,  development  of  explicit  formulae  for  ^(^)  limited  to 

certain  queueing  models  and  costs  function.  We  shall  demonstrate  the 
procedure  by  deriving  the  formulae  for  the  M/M/1  model  and  linear  costs. 
The  procedure  is  applicable  however  to  other  cases. 

Consider  the  M/M/1  model,  specified  by  the  set  of  parameters 

the  service  rate  function  p(x)  = max{xp,p]  and  cost  func- 
tion 


(3.8) 


m(x)  = — [ s + qx] 
P 


Let 


(3.9) 


+ n + p - ((\^  + + n + p)^  - 

1 ? 1 /p 

+ Q'i  + P + P + ((>^^  + ct^  + p + p)  - Ux^p)  ' ] 

1 

j=0,l,2,.... 


The  following  procedure  is  limited  to  cases  where  4 ^k 
i ^ k.  This  limitation  presents  some  complications  in  development.  Other 
assumptions  may  lead  to  more  complicated  results. 


Suppose  we  have  for  a given  k,  0 < k < n 


(3.10) 


M^(k)  = 


s + q(x  + 


+ i;  c 


I, 


p j=k 


for  appropriate  sequences  of  constants  k - 0,...,n;  j 

Then,  following  (3.7) 


Q. 


(3.11) 

Furthermore,  according  to  (2.14) 


n 


+ SO, 


p j=k 


a. 


(3.12) 

where 


,,\-l'^P  Ir  , , , “k-1 

U — Ls  + q(x  + 

P a^.l+P 


= - 
P P 


s + q(x  + 


(3.13) 


°k-l  \ ^ P^-1 


■k-1  + P 

From  this  result  and  (2.12)  we  obtain  that 


(3.14) 


u 


_x+l 

Q:.  i+P  Ct,  , n %. 

k-1  _krl_  . £ S C.  ^ 


"^k-l-^P  P i=k  ^-?i 

/^x+1  _x+l 

_ a s c • ( — . JEaU 


where 


(3.15) 


k-1  -J 


'k-1,1  - - 5^.i)(Vi-5))  "J 


O],,  > j - k, . 


Moreover,  since 


(3.16)  ^k-l<?j-5k-l)«k-l-Ej>  ■ \ 


we  obtain 


a. 


(3.1T) 


k-1 


-k-l,j  - (\,  ,-x..)(i-5.)  +a  -a.  "k,j 


C,.  A > J 


r 


Finally,  introduce  the  constants 


(3.18) 


II 

c-l,k-l  ' ^ ^k-l,i 


and  combine  the  partial  results  obtained  above  then 


(3.19) 


'*’P 

Mk-l(x)=U  “k-l^x)  = 

1 f \-l"'^  " 

= — s + q(x  + + S C,  . . — ■ ■ 

P \ P j=k-l 


Suppose  we  assume  that  ® ® given  n,  which  means  that 


A(t)  = for  all  < t then 


(3.20) 


/ \ -H  5 

^n^x)  = - S + q(x  + + — 

n p \ p l-F 


Commencing  with  q = X and  C = 1 we  apply  (3.13)>  (3.1?)  ^nd 
n n n,n 

(3.18)  to  the  recursive  computation  of  the  coefficient  q.,  C , 

K K ^ X 

i = k,...,n,  for  all  0 < k < n. 

Suppose  now  that  we  have  1°*’  some  pair  j,k  j > k. 

Formula  (2.12)  cannot  be  used  then  for  the  derivation  of  (3.1^).  Instead 
of  (2.12)  we  have  to  apply  (2.13)  and  introduce  terms  of  the  form  x|  . 

The  introduction  of  such  terms  complicates  the  equations.  Further  develop- 
ment may  be  possible,  but  it  requires  the  use  of  expression  (2.9).  The 
procedure  described  here  may  be  proved  unpractical  in  cases  where  the  s 
are  not  distinct. 

In  some  cases  there  is  certain  freedom  in  the  selection  of  the  sequences 
gi  and  \.  This  is  the  case  where  fA(t),  t > 0}  is  designed  to  approximate 
some  given  deterministic  function  \(t).  In  such  cases,  whenever  it  is 


necessary  to  assume  that  for  1 ^ k it  Is  still  possible  to 
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